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BSTRACT  rCanllnua  an  racarca  afSa  II  naeaaaaqp  anU  UmMIr  hr  Mack  wikaaj 

Before  satellite  images  obtained,  with  multiple  image  sensors  can  be  used 
in  image  analysis,  corrections  must  be  Introduced  for  the  differences  in 
transfer  functions  of  these  sensors.  Methods  are  here  presented  for 
obtaining  the  required  Information  directly  from  the  statistics  of  the 
sensor  outputs.  The  assumption  is  made  that  the  probability  distribution 
of  the  scene  radiance  seen  by  ea  di  Image  sensor  Is  the  same.  Successful 
destrlping  of  LANDSAT  images  is  demonstrated. 
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ARTIFICIAL  INTELLIGENCE  LABORATORY 


DESTRIPING  SATELLITE  IMAGES 


B.  K.  P.  Horn  and  R.  J.  Woodham 


Abstract.  Before  satellite  Images  obtained  with  multiple  Image  sensors  can 
be  used  In  Image  analysis,  corrections  must  be  Introduced  for  the  differences 
In  transfer  functions  of  these  sensors.  Methods  are  here  presented  for  ob- 
taining the  required  Information  directly  from  the  statistics  of  the  sensor 
outputs.  The  assumption  Is  made  that  the  probability  distribution  of  the 
scene  radiance  seen  by  each  Image  sensor  Is  the  same.  Successful  destrlping 
of  LANDSAT  Images  Is  demonstrated. 


This  report  describes  research  done  at  the  Artificial  Intelligence  Labonatory 
of  the  Massachusetts  Institute  of  Technology.  Support  for  the  laboratory's 
artificial  Intelligence  research  Is  provided  In  part  by  the  Advanced  Research 
Projects  Agency  of  the  Department  of  Defense  under  Office  of  Naval  Research 
contract  N000l4-75-C-0643. 


1_.  Destripinq  of  iwages  obUined  usino  ■ultiol# 


An  Image  sensing  device  using  a single  photoelectric  sensor  which  Is 
mechanically  scanned  across  the  scene  produces  outstanding  digitized  Images 
since  sensitivity,  resolution  and  transfer  functions  ore  the  same  for  all 
points  In  the  Image.  Unfortunately,  such  a device  Is  limited  In  speed  by  the 
mechanical  movement.  More  Importantly,  It  Is  limited  In  speed  by  the  fact 
that  an  accurate  measurement  of  scene  radiance  requires  the  collection  of  an 
adequate  number  of  photons.  This  explains  the  preponderance  of  linear  arrays 
of  sensors  and  area  sensors  such  as  vidicons  which  are  otherwise  deficient 
because  of  geometric  distortions,  non-unifom  response,  non-uniform  resolution 
and  so  on. 

A compromise  can  be  struck,  where  a small  set  of  sensors  Is  mechanically 
scanned  to  collect  the  Image.  In  the  system  used  aboard  LAMDSAT,  for  example, 
each  spectral  band  Is  scanned  using  six  sensors  at  the  same  time.  Thus,  six 
lines  of  the  Image  are  produced  during  a single  sweep  of  the  mirror.  On  the 
next  sweep  the  satellite  has  advanced  In  Its  orbit  by  an  amount  which  allows 
the  same  set  of  sensors  to  pick  up  the  next  six  lines  of  the  Image  CH. 

Unfortunately,  the  sensors  do  not  have  Identical  transfer  functions.  As 
a result.  Images  produced  In  this  fashion  show  undesirable,  regular  "striping”, 
This  effect  can  be  resNived  If  the  transfer  functions  are  accurately  known, 
since  one  could  then  compute  scene  radiance  from  the  sensor  output  using 
the  Inverses  of  these  transfer  functions.  The  sensors  used  In  older  equipment 
In  particular  have  time-varying  behavior.  Photomultipliers,  for  example,  show 
a drift  In  both  gain  and  offset  (dark  current)  due  to  small  changes  In  the 
material  of  the  dynodes  used  In  the  electron  multiplier  stages  and  temperature 
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varlatlons. 

If  a rafaranca  objact  containing  all  scana  radlancas  of  Intarast  wera  | 

In  the  scana,  ona  could  racallbrata  tha  sansors  continuously.  This  Is  diffi- 
cult to  arranga.  An  altamativa  Is  tha  scanning  of  a grio^  wedge  placed  over 
a light  source  at  tha  and  of  every  scan  line.  This,  In  fact.  Is  what  is  dona  I 

aboard  LANOSAT.  Tha  results  are  used  to  astinata  tha  gains  and  offsets  of  the  | 

1 

sensors.  Tha  digital  data  produced  froa  tha  raw  satellite  signals  Is  corrected  j 

using  this  Infonaatlon  [1].  | 

Unfortunately,  ona  finds  that  tha  striping  effect  Is  not  rewovad  In  this  j 

I 

fashion;  tha  reasons  for  this  are  not  entirely  clear.  One  cause  appears  to  | 

be  the  use  of  the  calibration  data  as  a means  of  adjusting  gain  and  offset 
so  that  each  sensor  Is  related  to  Its  prafllght  condition.  Slight  changes  In 

) 5 

the  light  source,  the  gray  wedge  and  the  geometry  of  Imaging  Introduce  drifts 
which  are  not  compensated  for.  Another  reason  Is  related  to  the  fact  that 
photomultipliers  are  somewhat  nonlinear  and  have  a response  which  depends  on 
their  exposure  history.  Modem  devices  using  solid  state  photodiodes  do  not  { 

suffer  from  these  problems. 

The  methods  explored  here  for  destrlping  Images  are  based  on  the  assump- 
tion that  each  sensor  Is  exposed  to  scene  radiances  with  approximately  the  ^ 
same  probability  distribution.  The  sensor  values  can  then  be  modified  so  that 
each  one  Is  related  In  the  same  way  to  the  actual  scene  radiance.  The  in- 
formation required  to  perfoim  this  modification  Is  extracted  from  statistics 


I 


2. 


A $1«p1e  Method  for  IlMar  trani 


If  the  Image  sensors  are  linear  and  time  Invariant,  a simple  method  can 
be  used  to  reduce  striping.  The  sensor  output,  x',  can  be  written  as  a func- 
tion of  the  scene  radiance,  x,  as  follows: 

x'  ■ f(x)  ■ a + b • X 

1 

j 

Each  sensor  has  Its  own,  fixed  values  of  offset,  a,  and  gain,  b.  If  these 

are  known,  the  scene  radiance  can  be  calculated  using  the  Inverse  of  the  I 

transfer  function. 


X • g(x')  ■ (x*  - a)/b 


If  this  Is  done  for  each  sensor  In  turn,  striping  effects  will  be  removed. 

The  required  constants  for  each  sensor  can  be  determined  If  a calibration 
object  containing  two  or  more  known  scene  radiance  values  Is  available  In  the 
scanned  scene.  If  such  a calibration  object  Is  not  available  one  can  estimate 
the  (relative)  values  of  gain  and  offset  using  simple  statistics  of  observed 
sensor  values.  Each  sensor  sees  a subimage  consisting  of  every  nth  line 
(when  n sensors  are  used).  The  complete  Image  Is  formed  by  Interlacing  these 
subimages.  It  seems  reasonable  to  suppose  that,  for  a large  enough  Image, 
each  subimage  has  approximately  the  same  probability  distribution  of  scene 
radiance  values.  One  would  not  expect  a particular  subimage  to  contain  many 
more  values  In  a particular  range  of  scene  radiances  than  another  subimage. 

If  this  assumption  Is  correct,  then  the  gain  and  offset  constants  can  be 
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estlmated  from  the  mean  and  standard  deviation  of  the  measured  sensor  output 
values.  If  the  mean  of  the  scene  radiance  Is  m and  the  standard  deviation  Is 
o,  then  the  mean  of  the  sensor  output  will  be  u ' ■ a b u and  the  standard 
deviation  of  the  sensor  output  a*  ■ b o.  Then, 

b » a'/a 

and 


a ■ (y'  o - p o' )/o 

Clearly,  It  Is  not  reasonable  to  assume  that  one  can  find  the  absolute  values 
of  the  mean  and  standard  deviation  of  the  actual  scene  radiance.  Fortunately, 
for  destrlping  purposes  only  relative  values  are  Important.  That  Is,  one  can 
use  the  mean  and  standard  deviation  of  the  sensor  outputs  for  the  whole  Image 
In  place  of  the  mean  and  standard  deviation  of  the  scene  radiance.  Naturally 
now  the  results  will  not  be  scene  radiance  values.  The  striping  however  will 
be  removed  since  each  subimage  now  has  the  same  mean  and  standard  deviation, 
and.  If  the  assumption  Introduced  earlier  applies,  the  same  linear  relation- 
ship to  scene  radiance. 

Note  that  one  can  relax  the  assumption  about  the  relationship  of  the 
subimages.  Here  It  Is  not  necessary  that  they  have  the  same  probability  dis- 
tribution of  scene  radiance,  only  that  their  means  and  standard  deviations  be 
the  same.  This  simple  method  has  been  applied  by  some  users  of  LANOSAT  daU 
[2.3]. 


3.  Shortcomings  of  the  simple  method 


We  have  found  this  method  to  be  only  partially  successful  In  destrlping 
LANDSAT  Images.  One  reason  for  this  may  be  that  out  of  a range  of  128  possible 
sensor  outputs  a range  of  only  around  30  values  correspond  to  noimal  scene 
radiance  values.  Low  values  are  not  found  In  short  wavelength  bands  because 
of  light  scatter  In  the  air.  Conversely,  large  values  correspond  to  cloud, 
snow  and  Ice,  and  scene  radiance  values  of  such  areas  often  exceed  the  highest 
available  sensor  output  values  and  so  result  In  clipping.  Clipping  of  sensor 
values  corresponding  to  low  scene  radiances  also  occurs  at  times  In  the  long 
wavelength  bands  due  to  negative  sensor  offsets.  Both  of  these  nonlinear  ef- 
fects will  Introduce  skew  Into  the  calculation  of  means  and  standard  deviations 


One  may  alleviate  this  problem  by  removing  sensor  values  outside  a certain 
range  from  consideration.  While  slightly  better  results  are  obtained  In  this 
fashion,  it  Is  dear  that  the  arbitrarily  selected  thresholds  needed  Introduce 
biases  of  their  own.  This  later  effect  can  be  dealt  with  by  eliminating  the 
same  fraction  of  sensor  values  from  the  low  end  of  the  output  of  each  sensor. 
Similarly,  a fixed  fraction  of  sensor  values  Is  removed  from  the  high  end. 

Even  with  this  refinement,  results  are  not  entirely  satisfactory.  Super- 
ficially, It  appears  that  different  gains  and  offsets  are  appropriate  for 
different  scene  radiance  ranges.  That  Is,  the  sensor  transfer  curves  are  some- 
what nonlinear.  We  thus  devised  a method  which  deals  with  this  problem  directly 
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4.  Preliminary  considerations. 

Consider  a random  variable  X with  probability  density  function  p(x). 
The  function  p(x}  Is  non-negative  and  satisfies 


f’  P(«) 


dx  » 1 


The  probability  density  function  p(x)  can  be  estimated  from  a large  number  N 
of  observations  of  the  random  variable  X.  If  n of  these  measurements  fall  In 
the  Interval  [x,  x + «x],  then  n/N  tends  to  p(x)  «x  as  N becomes  very  large 
and  6x  small  (In  a fashion  which  allows  N 6x»  and  thus  n,  to  become  large 
also). 

The  cumulative  probability  density  function  P(x)  Is  defined  as 


P(x)  - p(t) 


This  function  Is  monotonically  non-decreasing  since  p(x)  Is  non-negative.  P(x) 
represents  the  probability  that  the  random  variable  X takes  on  a value  less 
than  or  equal  to  x. 

Now  consider  observing  the  random  variable  X by  means  of  a transducer 
with  transfer  function  f(x).  Its  output  can  be  thought  of  as  a new  random 
variable  X',  say,  with  a probability  density  function  p'(x').  This  function 
Is  related  to  the  probability  density  function  p(x)  of  the  original  random 
variable  X,  In  a fashion  which  depends  on  the  transfer  function  x'  ■ f(x). 

It  Is  easiest  to  develop  this  relationship  In  terms  of  the  cumulative  distrl- 
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butlon  functions  P(x)  and  P'(x')  where 


I P'(x')  - / p'(t)  dt 

^ i ’ 

I I 

I If  x*  lies  In  a range  R*  when  x lies  In  the  range  R,  then  clearly. 


j p'(x')  dx'  » I p(x)  dx 
R‘  R 


Now  assume  that  f(x)  Is  monotonically  non-decreasing. 


X S Xq  Is  mapped  Into  the  range  x*  ^ ^(Xq). 


Then  the  range 


P‘Cf(x)]  - P(x) 

As  a result  one  can  determine  the  transfer  function  f(x)  If  the  cumulative 
probability  density  functions  P(x)  and  P'(x)  are  known  and  If  the  latter  has 
an  Inverse.  Then, 

f(x)  « (P*)“^  P(x) 

i 

1 

! 

If  P*  Is  monotonically  Increasing,  the  required  Inverse  will  exist.  Difficulties  • 

will  be  encountered  only  when  P'(x')  is  constant  over  a certain  range.  That  Is,  ' 

i 

If  P'(x')  ■ c [and  hence  p'(x')  ■ 0]  for  x'  e [xj,  x^].  Then,  If  P(x)  * c,  one  j 

can  say  only  that  f(x)  e [xl,  xi].  | 

There  are  two  possible  causes  of  this  problem.  First  It  may  be  that  f(x)  I 


actually  has  a discontinuity.  In  this  case,  one  correctly  finds  a jump  from 
Xi  to  X2  In  the  solution.  The  other  possibility  Is  more  serious.  If  p'(x')  » 
0 because  p(x)  » 0 [where  x'  * f(x)  as  before],  then  the  transfer  function 
f(x)  cannot  be  found  In  the  specified  range  because.  In  essence,  no  Inputs  are 
available  to  test  It  In  this  range.  The  Information  to  recover  f(x)  there 
Is  thus  not  available. 

Note,  however,  that  If  the  Inputs  to  the  transducer  are  In  fact  character 
Ized  by  the  given  probability  density  function,  then  our  lack  of  knowledge 
of  the  transfer  function  In  the  specified  range  Is  of  no  consequence  since 
no  inputs  fall  in  this  range  anyway. 

To  calculate  scene  radiance  from  sensor  values,  we  actually  need  the  In- 
verse g(x')  of  the  transfer  function.  This  can  be  found  just  as  easily.  If, 

P'(x')  - P [g(x')] 

Then 

g(x')  = P"^  P'{x') 

The  same  considerations  regarding  the  existence  of  the  Inverse  P**^  apply  here 
as  those  discussed  regarding  the  existence  of  the  Inverse  (P')~^  • All  these 
special  case  problems  are  avoided  If  the  cumulative  probability  distribution 
functions  are  monotonically  Increasing. 

The  method  shown  here  for  finding  the  transfer  function  of  a transducer 
(or  Its  inverse)  is  based  on  the  same  analysis  as  that  used  to  design  a genera 
tor  of  pseudo-random  numbers  with  a desired  probability  distribution  function 


p'(x)  when  a generator  Is  available  which  produces  pseudo-random  numbers  with 
known  probability  distribution  function  p(x)  [4]. 

A graphical  Illustration  of  the  relationships  discussed  may  be  found  In 
Figure  1.  The  dotted  line  suggests  how  one  may  determine  the  transfer  func- 
tion value*  f(x),  given  a scene  radiance  value  x.  Conversely*  the  same  dotted 
line  may  be  followed  In  the  reverse  direction  to  find  the  value*  g(x')*  of  the 
Inverse  function  from  a given  value  of  the  transducer  output  x'. 
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5.  Estimation  from  a finite  number  of  samples. 

To  apply  this  method  to  determine  the  transfer  function  of  a real  trans- 
ducer, the  cumulative  probability  density  functions  must  be  determined  from  a 

model  of  the  underlying  process  generating  the  random  variables  or  estimated 

I > 

from  frequencies  of  observed  occurrence  using  a finite  number  of  samples.  In 

the  latter  case  an  uncertainty  will  be  found  In  the  estimation  of  the  proba- 
bilities which  will  be  Inversely  proportional  to  the  square  iDot  of  the  number 
of  samples  falling  In  a particular  Interval. 

In  particular,  the  sample  deviation  of  the  estimate  of  the  cumulative 
;probab111ty  P(x)  obtained  from  a total  of  N measurements  Is 

aP  = /Pd  - P)/N 

The  following  approximate  analysis  may  be  helpful:  Let  af  be  the  uncertainty 
In  estimating  the  transfer  function,  x'  = f(x),  resulting  from  the  uncertainty  aP 
In  estimating  P(x},  then, 

P'[f(x)  + af]  = P(x)  + oP 

If  P'(x')  Is  reasonably  well  behaved  then  the  left-hand  side  can  be  expanded, 
P'[f(x)]  + p'(x')  af  z P(x)  + aP 

using  the  fact  that  p'  Is  the  first  derivative  of  P*.  Clearly  then 


p'(x')  of  = oP 


The  uncertainty  In  estimating  the  transfer  function  x'  ■ f(x)  Is  thus  highest 
Mhere  p'(x')  Is  small.  In  fact,  we  have  already  shown  that  f(x)  cannot  be 
uniquely  determined  when  p'(x')  > 0. 

One  can  similarly  show  that 


where  og  Is  the  uncertainty  In  estimating  the  Inverse  transfer  function, 

X ■ g{x'),  due  to  uncertainty  aP  In  estimating  P'(x').  Not  surprisingly,  un 
certainty  In  estimating  the  Inverse  transfer  function  Is  highest  where  p(x) 


The  probability  density  functions  p(x)  and  p'{x')  are  related  by 


when  the  Indicated  derivatives  exist 


-12- 


I i 


6.  Transducer  with  discrete  output  values. 

Essentially  the  same  method  may  be  used  If  the  transducer  produces  dis- 
crete outputs.  Consider,  for  example,  a case  where  the  Input  range  can  be 
broken  up  Into  a number  of  Intervals  such  that 


f(x)  » 1 If  X e [x^, 


The  probability  density  function  of  the  output  of  the  transducer  Is  then  dis- 
crete and. 


11m  /^+1  ■ ® 


'1 


Clearly,  p|  ^ 0 and 


I Pi-1 

i ■ -«B 


The  cumulative  probability  density  function  can  be  defined  as  follows 


If  f(x)  Is  monotonically  non-decreasing,  then  the  same  argument  applied  In 
the  continuous  case,  leads  again  to 
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If  P*  can  be  Inverted,  the  transfer  function  can  be  found  using 
f(x)  - (P')‘^  P(x) 

The  only  difference  Is  that  here  f(x)  Is  a function  from  a continuous  range 
to  a discrete  domain.  Naturally,  when  one  finds  the  Inverse  of  the  transfer 
function,  g(x'),  using  these  methods,  one  has  to  accept  the  fact  that  the 
actual  value  of  x cannot  be  recovered,  only  a range  [x^, 


A particular  sensor » In  a system  using  n sensors,  sees  a sublmage  contain- 
ing every  n^^  line.  To  apply  the  methods  developed  here  one  has  to  assume  that 
each  sensor  Is  exposed  to  scene  radiances  with  similar  probability  distribu- 
tions. If  the  Image  Is  large  enough.  It  is  very  unlikely  that  one  sensor  will 
see  substantially  fewer  or  more  scene  radiances  In  a particular  range.  A 
sensor's  properties  can  then  be  estimated  from  the  statistics  of  Its  outputs. 

Since  the  probability  distribution  function  of  the  actual  scene  radiance 
Is  not  available,  only  relative  adjustments  can  be  made.  That  Is,  the  proba- 
bility distribution  function  of  sensor  outputs  for  the  whole  Image  Is  used  as 
a reference.  Consequently,  application  of  the  Inverse  functions  so  determined 
to  each  subimage  will  not  produce  scene  radiances.  It  will,  however,  result 
In  an  Image  In  which  each  image  gray  level  Is  related  to  the  scene  radiance 
In  the  same  way.  Thus,  striping  will  have  been  removed. 

A graphical  Illustration  of  the  method  Is  shown  In  Figure  2.  Note  that 
here  both  P(x)  and  P'(x')  are  discrete.  This  results  In  a small  problem 
when  the  transfer  function  f(x)  Is  to  be  found  using  data  from  a real  Image. 

For  perfect  data,  there  always  Is  a value  x'  for  every  value  x,  such  that 
P'(x')  = P(x).  As  Indicated  by  the  dotted  line  In  the  figure,  this  may  not  be 
the  case  when  the  data  Is  obtained  from  a finite  number  of  samples  obtained 
by  different  sensors.  The  best  one  can  do  then  Is  to  find  a value  x'  such  that 

P'(x')  i P(x)  < P'{x'  + 1) 


Similarly,  In  determining  the  Inverse  transfer  function  g(x'),  one  can  do  no 
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better  than  finding  a value  x such  that 


i 


P(x)  S ?'(x')  < P(x  + 1) 


What  value  should  be  used  In  the  lookup  table  for  destrlping?  One  might  argue 

that  some  values  should  be  translated  to  x,  others  to  (x  1).  If  this  Is  done 
In  the  appropriate  fashion,  the  histogram  of  gray  levels  In  the  destrlped  Image 
will  equal  the  histogram  of  gray  levels  In  the  raw  Image.  It  Is  difficult 
however  to  decide  which  points  should  receive  one  value  and  which  the  other. 

If  the  selection  Is  based  on  a random  number  generator,  additional  noise  will 
be  Introduced.  In  any  case  a change  smaller  than  one  gray  level  Is  usually 
negligible  considering  the  limitations  of  radiometric  accuracy  In  the  Imaging 
system.  The  algorithm  we  employed,  arbitrarily  uses  the  smaller  value  of  the 
two  possible  ones.  Note  that  the  histogram  of  the  destrlped  Image  will  not 
be  exactly  equal  to  the  histogram  of  the  raw  Image  when  this  Is  done. 


i 


I 

1 
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8.  Details  of  the  algorithws. 

The  first  step  Is  the  detereil nation  of  a cumulative  histogram  of  sensor 
values  for  the  whole  Image  as  a reference.  Let  there  be  H(x)  occurrences  of 
sensor  outputs  less  than  or  equal  to  x out  of  a total  of  N values.  Now  for 
. the  subimage  produced  by  sensor  1,  one  calculates  a similar  cumulative  histo- 
gram. Let  H^(x')  be  the  number  of  sensor  outputs  less  than  or  equal  to  x', 
produced  by  sensor  1»  out  of  a total  of  N^  values.  Here. 

n 

N - I N. 

1-1  ’ 

where  n Is  the  number  of  sensors. 

A lookup  table  g(x')  Is  now  constructed  by  applying  the  Inverse  of  the 
function  H(x)  to  H^(x').  This  lookup  table  Is  then  used  to  modify  all  the 
sensor  values  produced  by  sensor  1.  The  Inverse  can  be  calculated  relatively 
easily  since  H(x)  Is  a monotonically  non-decreasing  function.  The  lookup 
table  value  g(x')  Is  the  smallest  number  x such  that 

N^  H(x)  ^ N H^(x')  < N^  H(x  + 1) 

This  process  Is  repeated  for  each  sensor  In  turn,  until  all  Image  values  have 
been  modified  by  the  lookup  table  appropriate  to  the  sensor  with  which  they 
were  measured. 
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9.  Results. 

The  Image  used  for  the  Illustrations  here  Is  a part  of  LANDSAT-1  Image 
E-1078-09555  taken  1972/0ctober/9.  The  part  used  contains  364  lines  each  of 
430  pixels.  The  worst  striping  In  this  case  occurred  In  band  6 (.7)i  to  .Su) 
and  band  7 (.8y  to  l.ly),  so  the  discussion  will  concentrate  on  these  two. 

The  data  transmitted  from  the  spacecraft  Is  pseudo-logarlthmically  compressed 
and  converted  to  six  bit  numbers.  On  the  ground  the  compression  Is  removed 
using  a lookup  table.  The  result  Is  a seven  bit  number.  Naturally,  about 
half  of  the  possible  seven  bit  numbers  never  occur  in  a given  subimage.  Fur- 
thermore, scene  radiance  corresponding  to  normal  surfaces  (Other  than  Ice, 
snow  or  cloud)  occupy  only  the  lower  25%  to  35%  of  the  available  range.  As  a 


result,  the  total  range  of  scene  radiances  of  Interest  corresponds  to  a very 
small  set  of  distinct  sensor  output  values.  This  contributes  a little  to  diffi- 
culties In  completely  removing  striping  effects. 

Comparison  of  Figure  3,  the  original  band  6 Image,  with  Figure  4,  the  pro- 
cessed version,  shows  that  the  heavy,  regular  striping  Is  removed  by  the  pro- 
cessing described  here.  It  Is  Instructive  to  Inspect  the  lookup  tables  used 
for  each  sensor.  These  are  shown  as  six  subfigures  In  Figure  5.  The  short 
horizontal  sections  In  the  transfer  function  correspond  to  sensor  values  which 
never  occur  as  a result  of  the  decompression  table  lookup.  The  Inverse  trans- 
fer functions  shown  In  Figure  5 appear  to  evidence  some  nonl Inearl  ties.  This 
explains  why  the  simple  destrlping  technique  described  earlier  does  not  do  as 
well. 

Comparison  of  Figure  6,  the  original  band  7 Image,  with  Figure  7,  the 
processed  version,  shows  that  the  regular  striping  Is  removed  here  as  well. 


i 


I 


i 

! 


M 

The  very  light  striping  rcMlnlng  In  shadow  areas  amounts  to  fluctuations  of  only 
one  or  two  pixel  values,  but  may  be  discemable  because  of  the  human  observer's 
sensitivity  to  small  differences  of  reflectance  In  dark  areas.  Inspection  of  the 
histograms  of  raw  Image  data  reveals  uneven  Intervals  In  the  analog-to-digital 
.convertors  used,  which  contribute  to  this  effect  (In  fact  some  codes  near  major 
transitions  never  occur).  The  relative  smoothness  of  the  lookup  tables  used  In 
destrlping  this  band,  shown  In  Figure  8,  Is  In  part  due  to  the  fact  that  no 
compression  and  decompression  Is  performed  on  data  for  this  band;  six  bit 
numbers  linearly  related  to  sensor  output  being  transmitted  directly.  The 
Inverse  transfer  functions  shown  In  Figure  8 appear  to  be  fairly  linear,  which 
Is  probably  a result  of  the  linearity  of  the  silicon  photodiodes  used  for  this 
\ Infrared  band.  One  would  expect  then  that  the  simple  destrlping  method  would 

I 

be  fairly  successful  for  this  band. 

We  experienced  some  difficulties  due  to  the  normalization  processing  per- 
formed on  the  raw  LANDSAT  data.  It  may  be  useful  to  provide  users  of  LANDSAT 
tapes  optionally  with  the  original  data.  Preliminary  results  Indicate  that 
slightly  better  destrlping  may  be  possible  using  the  raw  Image  sensor  values. 

It  Is  also  unfortunate  that  areas  of  high  reflectance  produce  scene  radiances 


which  saturate  the  Imaging  system.  If  this  was  not  the  case.  Image  sensor 
outputs  corresponding  to  areas  covered  by  thick  clouds  could  be  used  In  call- 
bratlon  of  absolute  reflectances  as  well  as  In  normalization  for  destrlping 
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10.  Relation  to  Histogram  Nomallzatlon  Methods. 

A number  of  techniques  aimed  at  the  enhancement  of  Images  Intended  for 
human  viewing  are  based  on  manipulation  of  the  gray  level  histogram.  Some*  for 
example,  transform  the  histogram  Into  one  considered  more  desirable  [5,  6,  7, 

8,  9],  either  flat  or  "hyperbolic".  The  suggestion  has  been  made  that  such 
techniques  may  allow  for  changes  In  sensor  characteristics  or  scene  Illumination 
[6,  8,  10].  Several  suggestions  have  been  made  regarding  the  difficulty  In- 
dicated In  Figure  2 relating  to  the  mismatch  of  two  cumulative  histograms. 

One  can,  for  example.  Increase  the  apparent  fineness  of  quantization  by  taking 
Into  account  the  context  of  each  picture  cell  [7].  If  each  picture  cell  has  a 
maximum  gray  level  value  m and  k neighbors,  one  may  multiply  Its  gray  level 
value  by  mk  and  subtract  (or  add)  the  sum  of  the  neighboring  gray  level  values. 
In  this  way  rank-ordering  of  picture  cell  values  Is  preserved,  while  the  number 
of  possible  gray  level  values  Is  vastly  Increased,  bringing  the  distributions 
closer  to  those  found  In  the  continuous  case. 

Recent  work  on  matching  of  Images  obtained  by  one  sensor  at  different 
times  Is  perhaps  most  closely  related  to  our  work  here  on  destrlping  [10].  This 
sort  of  approach  has  however  not  yet  found  much  of  a following  In  the  remote 
sensing  comnunlty  [11].  We  believe  our  application  of  histogram  "equalization" 
to  subimages  obtained  using  sensors  of  a multiple-sensor  system  Is  novel. 
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11.  Conclusions. 

The  folloiring  can  be  acconpllshed  using  the  method  presented  here: 

1.  The  effect  of  differences  In  the  transfer  functions  of  the 
different  sensors  Is  removed.  All  gray  levels  are  then 
related  In  the  same  *iay  to  the  original  scene  radiance  values. 

2.  The  overall  tone>sca1e  Is  preserved.  That  Is*  the  histogram 
of  gray  levels  of  the  destrlped  Image  Is  (approximately) 

the  same  as  the  histogram  of  the  row  Image.  There  Is  no  loss 
In  resolution*  nor  Is  the  noise  level  Increased. 

3.  The  computation  requires  only  two  serial  passes  over  the  whole 
Image:  one  to  collect  the  relevant  histograms*  another  to 
apply  the  Inverse  transducer  function  represented  as  a simple 
lookup  table. 


J 
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Figure  1. 

, Figure  2. 

Figure  3. 

Figure  4. 
Figure  5. 

Figure  6. 

Figure  7. 
Figure  8. 


The  transfer  function,  x*  ■ f(x),  of  a transducer  can  be  found  If 
the  cumulative  probability  distribution  functions  of  Its  Input,  x, 
and  Its  output,  x',  are  known. 

When  both  cumulative  probability  distribution  functions  are  discrete, 
there  may  not  be  a value  of  x'  so  that  P(x')  ■ P(x}.  Some  rule  must 
be  adopted  to  deal  with  this. 

Original,  unrectified  Image  of  band  6 (.7vi  to  .8y)  output  of  LANDSAT. 
Notice  the  heavy,  regular  striping. 

Destrlped  Image  of  band  6 output. 

The  destrlping  lookup  tables  fior  band  6 — Inverse  transfer  functions 
for  the  six  Image  sensors.  The  noni Inearl ties  of  the  transducers  are 
apparent. 

Original,  unrectified  Image  of  band  7 (.8ti  to  l.ly)  output  of  LANDSAT. 
Notice  regular  striping. 

Destrlped  Image  of  band  7 output. 

The  destrlping  lookup  tables  for  band  7 — Inverse  transfer  functions 
for  the  six  Image  sensors.  These  transducers  appear  to  be  fairly 
linear,  differing  mostly  only  as  regards  gain  or  amplification. 


FIGURE  3 


FIGURE  6 


